﻿using innovations.util.exts.mathdotnet;
using innovations.util.exts.msft;
using MathNet.Numerics.LinearAlgebra.Double;
using MathNet.Numerics.LinearAlgebra.Generic;

namespace innovations.ml.core.models
{
    public class LinearRegressionModel : Model, IModel
    {
        public LinearRegressionModel(Matrix<double> x, Vector<double> y, bool includesOnesColumn) : base(x, y, includesOnesColumn) { }
        public virtual void ComputeCost(bool includeGradient = false)
        {
            double m = Y.Count; // number of training examples
            Vector<double> sigmoidValue = Sigmoid.Compute((X).Multiply(Theta));
            int r = 1;
            J = (1 / m) * (-Y.PointwiseMultiply(sigmoidValue.PointwiseLog()) - (r.PointwiseSubtraction(Y).PointwiseMultiply((r.PointwiseSubtraction(sigmoidValue)).PointwiseLog()))).Sum();
            if (includeGradient)
                ComputeGradient();
        }

        protected virtual void ComputeGradient()
        {
            double m = Y.Count;
            Vector<double> temp = new DenseVector(X.ColumnCount);
            for (int j = 0; j < X.ColumnCount; j++)
            {
                Vector<double> vector = ((X).Multiply(Theta)).Subtract(Y);
                Vector<double> xVector = X.Column(j);
                temp[j] = Theta[j] - (Solver.Alpha * (1 / m) * (vector.DotProduct(xVector)));
            }
            Theta = temp;
        }

        /// <summary>
        /// Used by ALGLIB.  All parameters are ignored.
        /// </summary>
        /// <param name="theta">ignored</param>
        /// <param name="j">ignored</param>
        /// <param name="gradient">ignored</param>
        /// <param name="obj">ignored</param>
        public virtual void ComputeCostAndGradient(double[] theta, ref double j, double[] gradient, object obj)
        {
            Theta = new DenseVector(theta);
            ComputeCost(true);
            j = J;
            // dm.Gradient = gradient doesn't work with L_BFGS solver in ALGLIB.
            // Therefore, we have to iterate through eack gradient element 
            // and assign to dm.Gradient.
            for (int i = 0; i < Gradient.Count; i++)
                gradient[i] = Gradient[i];
        }
    }
}